function h2dinit(params)

  % ****************************************************************
  %  output info
  % ****************************************************************

  % summary
  disp(sprintf('Grid Size   : (%d, %d)', params.Nx, params.Ny));
  disp(sprintf('Grid Length : %f', params.Lx));
  disp(sprintf('Cell Length : %f', params.Dx));
  disp(sprintf('Target Time : %f', params.Tfinal));
  disp(sprintf('Time steps  : %d',  params.Tsteps));
  disp(sprintf('Hydro Size  : %6.2fMB', 3*params.Nx.*params.Ny.*5.*params.Tsteps.*64/8/1024/1024));



  % ****************************************************************
  %  prepare plotting
  % ****************************************************************

  % prepare figures for plotting
  close all;

  % load color maps
  load('cm', 'mycmap');
  global cmSize;
  cmSize = 512;

  % figures
  if (params.plotHydro)
    fig  = figure; 
    set(fig, 'Colormap', mycmap);
    set(fig, 'WindowStyle', 'docked');
  end
  refresh();
  pause(0.001);
  



  


  % ****************************************************************
  %  compute derivative kernels
  % ****************************************************************

  tic;

  Nx = params.Nx;
  Ny = params.Ny;
  Dx = params.Dx;
  Dy = params.Dy;

  % can upgrade on this -- trefethen pg 6

  % the derivative code assumes Nx = Ny so far

  global d1x d1y d2x d2y;
  global gd1x gd1y gd2x gd2y;

  d1x = zeros(Nx, Nx);
  d1y = zeros(Ny, Ny);
  d2x = zeros(Nx, Nx);
  d2y = zeros(Ny, Ny);

  % x derivatives
  for x = 1:Nx
    for y=1:Nx
      s1 = 0;
      s2 = 0;
      for k = -Nx/2 : Nx/2
        s1 = s1 + 2. * pi * i / Nx / Nx * k * exp(2 * pi * i * k * (x - y) / Nx);
        %s2 = s2 + (2. * pi * i / Nx * k)^2 / Nx * exp(2 * pi * i * k * (x - y) / Nx);
      end
      d1x(x,y) = real(s1)/Dx;
      %d2x(x,y) = real(s2)/Dx/Dx;
    end
  end

  % y derivatives
  for x = 1:Ny
    for y=1:Ny
      s1 = 0;
      s2 = 0;
      for k = -Ny/2 : Ny/2
        s1 = s1 + 2. * pi * i / Ny / Ny * k * exp(2 * pi * i * k * (x - y) / Ny);
        %s2 = s2 + (2. * pi * i / Ny * k)^2 / Ny * exp(2 * pi * i * k * (x - y) / Ny);
      end
      d1y(x,y) = real(s1)/Dy;
      %d2y(x,y) = real(s2)/Dy/Dy;
    end
  end


  tder = toc;
  disp(sprintf('Initializing derivative kernels, %.2f sec', toc));
end
